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1.  INTRODUCTION 

Defense  Meteorology  Satellite  Program  (DMSP)  Block  5D  satellites  have 
included  Topside  Ionospheric  Plasma  Monitor  (SSIE)  sensors  since  mid  1977^’^. 
Ionospheric  electron  density,  N  (840),  and  plasma  scale  height,  H  (840),  at 

“  P 

the  satellite  altitude,  nominally  840  km,  from  the  SSIE  sensors  are  routinely 

processed  at  the  Air  Force  Global  Weather  Center  (AFGWC)  for  specification  of 

the  state  of  the  earth's  ionosphere.  Figure  1  is  a  simplified  flow  diagram  of 

data  from  the  SSIE  sensors  within  AFGWC.  SSIE  data  are  stripped  from  the  DMSP 

satellite  telemetry  stream  and  written  to  a  temporary  file  (IEPREPFILE)  by  the 

AFGWC  DMSP  processing  system.  Program  SSIE  reads  the  raw  SSIE  data  from  this 

3  4 

file,  calculates  electron  densities  and  plasma  scale  heights  ’  ,  and  writes 

out  a  data  record  for  every  64  seconds  containing  an  average  Ng(840)  value  and 

Hp(840)  value  to  a  temporary  file  (IFLUXFILE).  Program  LDIPIE  reads  the  data 

records  from  this  file,  reformats  them,  and  writes  them  into  the 

c 

Astrogeophysical  Data  Base  (AGOB)  .  Program  SSIELD  reads  SSIE  Ng(840)  data 

for  the  most  recent  24  hours  from  the  AGDB,  constructs  10°  latitude  by  15° 

longitude  northern-hemisphere  grids  of  N  (840)  for  each  of  the  24  hours,  and 

writes  them  to  a  gridded-data  file  (FOURDGRIDS)  .  These  data  are  merged  with 

data  grids  from  the  f  F2  and  TEC  preprocessors  and  written  to  temporary  file 

o  2 

18  by  program  GRDOUT  as  observation  data  sets  .  The  Air  Weather  Service  (AWS) 

7  8 

4D  Ionospheric  Analysis  Model  ’  inputs  the  observation  data  sets  from 

temporary  file  18  and  produces  a  consistent  analysis  of  the  ionospheric 

electron  density  distribution  over  the  entire  northern  hemisphere. 

The  primary  objective  of  this  project  is  to  study  methods  for  improving 
the  use  of  the  SSIE  Ng(840)  and  Hp(840)  data  by  the  4D  model  system, 
identified  in  Figure  1  by  those  files  and  programs  within  the  dashed  box.  The 
studies  have  focused  on  possible  improvements  to  the  4D  model,  particularly  to 
its  internal  electron-density  profile  model,  and  to  the  various  4D  model 
system  preprocessors. 
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SSIE  telemetry  data 


Ng(840)  averages,  Hp(840) 


Data-type  18  records 


Ng(840)  grids 


Observation  data  sets 


2.  PROJECT  DEFINITION 


The  first  step  in  the  project  was  to  define  specific  tasks  to  pursue 
within  the  areas  of  research  defined  in  the  project  proposal  and  statement  of 
work.  This  study  resulted  in  the  following  list  of  tasks  to  be  completed 
during  the  project: 

1.  Determine  the  modifications  required  to  configure  the  4D  iono¬ 
spheric  analysis  model  for  a  full  global  analysis. 

2.  Investigate  the  use  of  the  profile  parametrization  used  in  the 

g 

IRI79  reference  ionosphere  as  a  replacement  for  the  Damon-Hartranft 
parametrization  currently  used  in  the  4D  model. 

3.  Derive  sets  of  height-profile  basis  functions  for  the  4D  model 
using  plasma  frequency  profiles  from  the  IR 179  reference  ionosphere. 

4.  Investigate  techniques  for  improving  the  preparation  of  data  grius 
of  Ng(840)  from  the  SSIE  Ne(840)  observations. 

5.  Work  with  personnel  from  AFGL  and  AFGWC  to  develop  improved  tech¬ 
niques  for  combining  TEC  and  SSIE  Ng(840)  data  into  a  single  specification  of 
the  topside  ionosphere. 

6.  Investigate  alternative  parametrizations  of  the  topside  ionosphere 
which  would  use  the  full  set  of  data  available  from  the  SSIE  sensors  and  all 
other  available  data,  including  Total  Electron  Content  (TEC)  data,  to  specify 
the  electron  density  profile  above  the  F2-layer  peak. 

7.  Begin  an  investigation  into  ways  to  use  data  from  the  SSIE  sensors, 
in  conjunction  with  other  observations,  for  modeling  the  ionospheric  sub- 
auroral  trough  as  part  of  the  4D  model  preprocessor  system. 

Tasks  1,  2,  and  3  have  been  completed,  and  work  has  begun  on  Tasks  4,  5, 
and  6.  These  tasks  fall  generally  into  three  main  studies:  1)  Topside 
Profile  Model  Study  (Tasks  2  and  6),  2)  Ionospheric  Data  Preprocessor  Improve¬ 
ments  Study  (Tasks4,  5,  and  7),  and  3)  4D  Model  Improvements  Study  (Tasks  1  and  3) . 
Work  completed  this  project  year  will  be  presented  separately  for  each  of 
these  three  study  areas. 
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3.  TOPSIDE  PROFILE  MODEL  STUDY 


This  study  includes  Tasks  2  and  4  of  the  project.  Tne  overall  objective 
of  tnis  study  is  to  investigate  improved  parametrizations  of  the  topside 
ionospheric  electron  density  profile  for  use  within  the  4D  model. 


3.1  Task  2:  IRI79  Topside  Parametrizatior, 


The  4D  uses  a  mode  1  of  the  electron  density  profile  in  tne  topsice 

ionosphere  to  interpolate  between  observed  values  of  f„F2  ana  h  F2  and  in-situ 

o  in 

observations  of  Ng(h)  and  Hp(h),  the  electron  density  and  plasma  (pressure) 
scale  height,  at  DMSP  altitudes.  The  profile  model  must,  therefore,  be 
designed  such  that  its  parameters  can  be  easily  adjusted  to  fit  the  input  aata 
and  behave  reasonably  from  h(nF2  to  2000  km.  In  this  context,  reasonable 
behavior  is  defined  as  follows: 


(1)  The  electron  density  scale  height,  defined  as 

1  -1 


H„(h)  -  - 


dS  Ne<"> 


!  I  ' 
v  1  I 


is  positive  throughout  the  topside. 

(2)  There  are  no  sharp  discontinuities  in  either  N  (h)  or  Hn(h)  (impor¬ 
tant  for  ray  traceing  applications). 


(3)  The  scale  height,  H  (h),  varies  smoothly  from  ■->  at  n_F2  through  the 
value  calculated  at  the  DMSP  altitude  (nominally  840  km)  to  values  at  2000  km 
consistent  with  the  expected  temperature  and  mean  ionic  mass  at  that  altitude. 

Currently,  the  4D  model  uses  the  modified  Damon-Hartranft  model 
denoted  the  DHR  model  (also  referred  to  as  the  R3TEC  model),  to  interpolate  in 
the  topside.  The  two  major  problems  with  this  model,  in  this  application,  are: 
i)  the  profile  assumes  that  the  only  ionic  constituent  at  all  altitudes  is  G+, 
and  ?.)  only  tnree  of  the  four  data  values  to  be  input  ( f Q F 2 ,  M3G00,  Ne(840), 
and  Hn ( 840 ) ;  can  be  fit  at  one  time.  An  interim  solution  to  the  second 
problem  was  developed  during  earlier  work  on  this  project  ,  but  the  resulting 

profiles  can  sometimes  exhibit  unreasonable  behavior  in  H  .'hi. 

n '  • 
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The  major  thrust  of  Task  2  was  to  investigate  the  IRI79  model  to  deter¬ 
mine  if  1)  the  IRI79  topside  parametrization  could  be  modified  for  use  within 
the  4D  model,  and  2)  the  profile  shape  provided  by  the  IRI79  model  is 
sufficiently  better  than  the  DHR  model  to  justify  replacement  of  the  DHR  model 
by  the  IRI79  model  in  the  4D. 

The  first  step  of  this  task  was  to  modify  the  IRI79  topside  profile 

parametrization  to  allow  the  profile  to  be  fit  to  input  topside  observations 

of  NJh)  and  H  (h).  The  IRI79  model  uses  an  exponential  function  for  the 
e  n  g 

topside  Ne(h)  profile  of  the  form 

N  (h)  =  N  F2  e-Y^  (2) 

e  m 


where 


Y(h)  =  R  [6nE1(h,394.5,  6  )  +  100  z, E^h^OO,  100)  ] 


(h-hmF2k 


(3) 


R  = 


(  1000-hmF2  ^ 

\  700  /  ’ 


Ej(h,a,b)  =  In 


1  *  w  (ji^) 

1  ♦  exp  (^) 


(h-h  F2) 

X(h)  =  - g -  +  300  -  6  , 


and  3,  n ,  z,  and  6  are  functions  of  various  geophysical  parameters.  The 
scale  height,  defined  as  in  Equation  1,  can  be  calculated  from  Equations  2  and 
3  to  give 


yh)  =  nE2(h,394.5,  3  )+  y  E£(h, 300,100)  -  Ij 


where 


E2(h,a,b)  =  R  ^  E^h.a.b)  = 


1  +  exp 


Ll  D 


An  iterative  technique  was  developed  to  fit  this  profile  to  input  values  o 
fQF2,  h^F2,  Ne(840),  and  y840)  by  adjusting  the  two  parameters,  •  and  c 
through 

=  |  [E2(h,  394.5,  3)  +  ±  E2 ( h , 300, 100) -l]  J  Hn(840)  1  (E 


where 


Y ( 840)  = 


-  Sn'E,  (h, 394.5, 


f -  lOOEjfh, 300,100) 


Ne(840) 


This  technique  was  successful  in  calculating  a  profile  which  reproduced  the 
'nput  data,  but  the  resulting  profiles  were  quite  often  unreasonable  away  from 
the  data  points,  at  times  producing  secondary  peaks  in  N  (h)  just  above  the  F2 
peak.  A  second  technique  was  then  developed  in  which  the  parameter  ,  S,  was 
adjusted  through 


i'(840)  +  (h-hmF2),-.H  -  100  E,  ;h,3G0, 100) 


h ,  394.5 


replacing  Equation  6,  but  this  technique  worked  no  better  than  the  first.  A 
study  of  the  IRI79  topside  was  then  begun  to  determine  the  source  of  the 
difficulties  in  producing  reasonable  profile  behavior  using  the  IRI79  topside 
parametrization. 

Since  the  desire  for  a  more  realistic  topside  profile  shape  was  the 

primary  reason  a  replacement  for  the  DHR  model  was  being  sought,  the  study 

focused  on  the  scale  height  profiles  produced  by  Equation  4.  The  Hn(h) 

profiles  thus  calculated  were  compared  with  the  results  of  an  analysis  of 

12 

Hn(h)  profiles  taken  from  Alouette  topside  soundings  .  Figure  2  shows  a 
representative  sample  of  daytime,  midlatitude  (dip  latitude  -45°),  solar 


O  200  400  600  800 

SCALE  HEIOHT  (KM) 

Figure  2.  Comparison  of  IRI79  and  Alouette  daytime,  midlatitude  scale  heights 

minimum  Hp(h)  profiles  from  the  IRI79  model  (solid  curve)  and  from  Alouette 
soundings  (dotted  curve  is  a  single  representative  ^n(h)  profile;  the  dashed 
curves  bracket  all  Hp(h)  profiles  for  the  Alouette  pass  from  which  the 
individual  profile  was  taken).  The  pattern  of  IRI79  Hp(h)  profile  departures 
from  Alouette  soundings  shown  in  Figure  2,  scale  heights  roughly  30 
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to  40%  too  low  just  above  the  F2  peak  and  up  to  300%  too  high  at  900  km,  was 
typical  for  all  daytime  comparisons.  Figure  3  shows  a  plot  of  scale  heights 
at  various  altitudes  as  a  function  of  dip  latitude  for  summer  day  and  summer 
night.  The  plot  for  daytime  shows  the  persistence  of  the  pattern  displayed  in 
Figure  2.  The  night  plot  shows  that  although  the  IRI79  scale  heights  are  not 
too  bad  at  lower  midlatitude  night,  the  model  departures  revert  to  their 
daytime  pattern  within  the  transition  into  the  light-ion  trough. 

A  final  comparison  was  made  between  TEC  values  calculated  by  Doth  the 
I R 1 79  and  DHR  models  and  observed  TEC  from  the  LlSAF  Sagamore  Hill  Observatory. 
The  model  profiles  were  calculated  using  ITS-78  model  values  for  fQF2  and 
M3000  modified  by  adding  the  difference  between  observed  and  ITS-78  values  at 
Wallops  Island,  which  is  near  the  Sagamore  Hill  ionospheric  penetration 
point.  Two  sets  of  data  were  used,  the  first  taken  from  monthly  median  data 
used  in  evaluation  of  the  Damon-Hartranft  model***  from  solar  cycle  20,  and 
five-day  mean  values  from  AF6WC  from  solar  cycle  21.  Table  1  shows  the 
results  of  this  comparison.  In  all  cases  the  DHR  model  was  appreciably  closer 
to  the  observed  TEC  than  the  IR 179  model.  While  not  as  extensive  a  comparison 
as  that  made  by  McNamara  and  Wilkinson*4,  the  data  used  are  fairly 
representative  of  median  midlatitude  behavior  and  favored  neither  model  a 
priori . 

TABLE  1.  Comparison  of  observed  and  model  TEC  (DHR  and  I R 1 79 ) . 


■j  r  O 

Total  Electron  Content  (xlOiD  el/m  ) 


Date 

SSN 

ut 

IT 

Observed 

DHR 

Model 

I R 1 79 

Model 

Jan  1968 

122 

0500 

0020 

6.6 

7.2 

(+11%) 

13.7 

(+108%) 

Jan  1968 

122 

2000 

1520 

41.0 

37.3 

Q0/ ' 

i +28%  j 

57.6 

(+  40% 

Mar  1968 

92 

1100 

0620 

6.4 

8.3 

12.7 

(+  95%) 

Mar  1968 

92 

1600 

1120 

33.0 

34.1 

i+  3%) 

44.5 

(+  35%) 

Jul  1969 

96 

0400 

2320 

12.0 

14.5 

(+21%) 

19.8 

(+  65%) 

Jul  1968 

96 

1800 

1320 

1S.0 

17.8 

{-  1%) 

21.4 

(+  19%) 

Sep  1968 

117 

0100 

2020 

14. C 

17.2 

(+23%) 

28.1 

(+101%) 

Seo  1968 

117 

1500 

1020 

24.0 

24.0 

(  0%) 

39.6 

(  +  65%) 
(  +  45%) 

Sep  1983 

86 

0000 

1920 

17.5 

19.7 

(+13%) 

25.4 

Sep  1983 

S6 

0600 

0120 

6.4 

5 . 6 

(-13%) 

9.5 

(+  48%) 

Sep  1983 

86 

1200 

0720 

11.0 

11.8 

(+  7%) 

17.7 

(+  61%) 

Sep  1983 

86 

1800 

1320 

17.2 

17.6 

(  +  2%) 

21.5 

(+  25%) 

In  conclusion,  it  was  decided  that  the  IRI79  model,  in  its  current  form, 
is  not  a  viable  replacement  for  the  DHR  model  for  the  following  reasons: 

1.  Although  the  topside  parametrization  used  in  the  IRI79  model  could 
be  modified  to  fit  input  observations  of  fQF2,  hmF2,  Ng(840),  and  8^(840),  the 
resulting  profiles  often  did  not  behave  reasonably,  occasionally  producing 
additional  Ne(h)  peaks  above  hmF2.  This  was  particuariy  true  if  the  input 
Ng(840)  or  Hn(840)  differed  significantly  from  that  calculated  from  IRI79 
using  only  fQF2  and  hmF2  as  inputs. 

2.  The  electron  density  scale  height,  Hn(h),  calculated  from  the  IRI79 
model  ranges  from  values  which  are  slightly  too  small  just  above  the  F 2  peak 
to  values  which  are  much  too  large  above  800  km,  particularly  during  the  day 
and  poleward  of  the  plasmapause  at  night. 

3.  The  top  of  the  I R 1 79  model  is  at  I0G0  kin,  and  the  40  model  requires 
a  profile  up  to  2000  km.  Extending  the  IR 179  profile  calculation  above  1000 
km  results  in  an  electron  distribution  with  an  effectively  infinite  scale 
height  above  roughly  1200  km. 

4.  The  sub-auroral  light-ion  trough,  a  major  feature  in  the  night irne 
topside  ionosphere,  is  not  present  in  the  IR 1 79  model. 

5.  Comparisons  of  TEC  calculated  from  profiles  generated  by  the  two 
models  showed  that,  at  best,  the  IRI79  model  performed  slightly  better  than 
the  DHR  model  only  on  occasion,  and  often  did  much  worse. 

3.2  Task  6:  Alternative  Topside  Parametrization 

Since  the  IR  1 79  model  proved  unsuitable  for  this  project,  work  was 
begun  on  Task  6,  developing  a  new  parametrization  for  the  topside  to  replace  the 
DHR  model.  The  basic  concept  that  evolved  is  to  develop  an  empirical  topside 
model  which  treats  the  ionosphere  profile  separately  above  and  below  the  height 
wnere  N(0  )=N(H  ),h-,  and  fits  the  two  orofiles  together  smoothly  at  n^  using  a 
spline  fit  procedure  which  produces  a  smooth  transition  in  Doth  N  (h)  and  H  (h) 
between  prof i le  sections.  This  requires  tnree  separate  models:  1)  an  Ng(h)  model 
for  the  height  region  h|nF2^h<hy,  2)  an  Ng(h)  model  for  hy<h^  2000,  and  3)  a 
model  for  the  benavior  of  h-.. 
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The  starting  point  for  the  profile  model  is  the  following: 


1)  h/2«h«hT 

Ne(h)  -  NfflF2exp  A  -  Zj  -  e  ’j 

vhere 

h-h  F2 
,  =  HL_ 

Z1  H1(h) 

Hl(h)  =  H0(  ^m’LT^  ln  8574 
<{>  :  magnetic  latitude 

LT  :  local  time,  and 

2)  hT<  h  ^2000  KM 

~z2 

Ng(h)  =  NTe 

where 

h-hj 

z2  =  T^TTTT 

h-h-r 

H2(h)  =  HTexp  -j~  . 

The  profile  in  the  lower  section  is  the  DHR  model,  which  is  a  reasonable  model 
for  a  predominantly  0  ionosphere.  The  upper  profile  is  an  exponential 
distribution  with  an  exponentially  varying  scale  height  with  several  free 
parameters,  Nj,  Hj,  and  B,  which  can  be  adjusted  as  necessary  to  fit  obser¬ 
vations. 

A  preliminary  model  has  been  developed  for  hj  based  on  analyses  of 
topside  sounder  data14  which  allows  hT  to  vary  as  a  function  of  magnetic 
latitude,  local  time,  season,  and  the  10.7  cm  solar  flux.  Figures  4a  and  4b 
show  examples  of  the  results  from  this  preliminary  model.  The  solid  lines  in 


the  plots  are  from  the  model,  and  the  dashed  lines  are  data  taken  from  Figures 

15 

9  and  11  in  Titheridge  (1976a).  Work  will  continue  on  this  model  during  the 

next  project  year.  The  solar  cycle  variation  of  hT  will  be  checked  against 

1 6 

recent  observations  ,  and  the  possibility  of  varying  the  location  of  the 

equatorward  edge  of  the  light-ion  trough  in  the  hT  model  as  a  function  of 
17 

magnetic  activity  or  from  observations  will  be  investigated. 
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Figure  4a.  Latitudinal  variation  of  hy  for  local  noon  and  midnight  (June). 

Solid  lines  are  from  topside  sounder  data^,  dotted  lines  a  ~e 
from  the  hT  model . 
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4.  IONOSPHERIC  DATA  PREPROCESSORS  STUDY 

Limited  investigations  into  analysis  stability  problems  with  the  4D 

model  led  to  the  development  of  front-end  data  preprocessors  for  f  F2  and  TEC 

18  19  20  0 

data  by  AFGWC/TSIS  ’  ’  .  These  preprocessors  use  simple  curve-fitting 

techniques  to  generate  equally  spaced  grids  of  f  F2  and  TEC  values  from  data 

0  21 

available  in  the  Astrogeophysical  Data  Base  (AGDB)  .  In  the  contract  prior 
to  the  current  one,  a  4D  model  preprocessor  system11  was  designed  around  these 
two  preprocessors,  programs  FORIER  (f  F2  data)  and  P0LY1  (TEC  data),  and  a 
third  preprocessor  for  SSIE  Ne(840)  data,  program  SSIELD6,  provided  by 
another  contractor. 

This  system,  tested  and  installed  at  AFGWC,  was  to  be  the  starting  point 
for  the  tasks  under  this  study  area.  Unfortunately,  as  of  the  start  of  this 
phase  of  the  project,  the  SSIELD  program  was  still  not  operational  at  AFGWC, 
and  was  neither  sufficiently  tested  nor  adequately  documented  to  be  of  any 
use.  Thus,  this  study  area  was  begun  by  assisting  AFGL/PHG  with  an  assessment 
of  the  SSIELD  program. 

4.1  Program  SSIELD  Assessment 

4.1.1  SSIELD  Overview 

At  the  request  of  AFGL/PHG,  tests  of  the  SSIELD  program  designed  to 
assess  the  program's  status  were  conducted  during  a  visit  to  AFGWC  in  late 
March,  1984. 

The  SSIELD  program  produces  equally  spaced  grids  of  N  (840)  from  SSIE 

fi 

data  in  the  AGDB  in  a  two-phase  process  .  First,  24  hours  of  SSIE  data  for  a 
single  DMSP  satellite  are  read  from  the  AGDB  and  are  used  to  construct 
equally  spaced  latitude/longitude  grids  for  the  ascending-  and  descending- 
node  local  times.  If  data  for  a  second  DMSP  satellite  are  available,  this 
step  is  repeated  for  data  from  the  second  satellite.  Second,  the  two  (or 
four)  local  time  grids  are  "spread"  to  other  local  times  by  a  weighted  inter¬ 
polation  scheme.  Output  from  the  program  consists  of  equally  spaced 
latitude/longitude  grids  tor  each  of  24  UT  hours. 

The  construction  of  the  local  time  grids,  or  maps,  is  accomplished  by  a 

22 

method  known  as  "krigina  ."  This  method  assumes  that  the  spatial  stochastic 
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field  described  by  the  data  can  be  described  by  the  sum  of  a  mean,  or  trend, 
field  and  a  fluctuation  field.  After  modeling  and  removing  the  trend,  tne 
fluctuation  field  is  modeled  in  terms  of  the  covariance  function  of  the 
process.  In  SSIELD,  the  trend  is  removed  by  fitting  the  data  for  a  given 
local  time  by  a  simple  polynomial  function  in  latitude  and  the  sine  of  the 
longitude^.  This  trend  is  removed  from  the  data,  and  the  resulting  fluctua¬ 
tion  field  is  used  to  estimate  the  covariance  function,  C.  In  estimating  C, 
two  major  assumptions  are  made: 

1.  The  fluctuation  field  is  homogeneous  and  isotropic.  Therefore  C  is 
a  function  only  of  the  separation  between  points,  and  not  a  function  of  point 
location  or  direction. 

2.  The  covariance  function  is  of  the  form 

C(d)  =  o2e''d  (8) 

2 

where  d  is  the  distance  between  two  points,  and  c  and  3  are  parameters  to 
estimate  from  the  fluctuation  field. 

2 

The  o  parameter  is  estimated  by  the  variance  in  the  fluctuation  field.  The  ■?. 
parameter  is  estimated  by  fitting  equation  8  to  a  discrete  covariance  function 
estimate  calculated  from  the  input  data.  Should  the  process  of  estimating 
fail,  the  program  uses  a  default  value  for  8  of  0.1.  The  trend  and  covariance 
function  estimates  are  then  used  to  produce  equally  spaced  local  time  maps. 

After  all  input  data  have  been  converted  to  local  time  maps,  a  model  of 
the  electron  density  at  840  km  is  used  to  interpolate  to  other  local  times. 
For  a  given  latitude/longitude  grid  point,  the  two  (four)  values  from  the 
local  time  maps,  weighted  by  error  estimates,  are  used  to  shift  the  local  time 
variation  curve  generated  from  the  Bent  model  for  that  grid  point.  This 
process  is  repeated  for  all  grid  points,  and  the  resulting  24  local  time  maps 
are  converted  to  24  UT  maps  and  written  out  to  file  FOURDGRIDS. 

The  electron  density  model  used  in  SSIELD  is  a  simple  spectral  expansion 

of  N  (840)  fields  extracted  from  electron  density  profiles  constructed  using 
e  23 

the  Bent  ionosphere  model  •  Four  sets  of  24  northern  hemisphere  Ng(840) 
grids,  one  grid  for  each  UT  hour  and  one  set  for  each  season,  were  fit  with 
associated  Legendre  function  polynomials  of  degree  9  and  order  2.  Interpola¬ 
tion  in  day-of-year  is  done  with  trigonometric  functions;  no  interpolation  is 
done  in  UT.  Although  this  model  represents  only  a  simple  snapshot  of  the  Bent 
model,  it  will  be  referred  to  as  the  "spectral  Bent  model"  in  the  remainder  of 
this  report. 
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4.1.2  Assessment  Test  Description 

A  series  of  tests  of  the  SSIELD  program  was  designed  to  assess  the 

program's  current  status  for  AF6L/PHG.  Nine  tests  were  planned  -  eight  using 

Np(840)  data  from  models  as  input  and  one  using  "live"  SSIE  data  from  the 

11 

AGDB.  Four  model  ionospheres  were  generated,  three  using  the  DHR  model 
(sunspot  numbers  10,  60,  and  120)  and  one  using  the  spectral  Bent  model  from 
SSIELD.  A  date  of  15  June  was  used  for  both  models,  to  avoid  a  negative 
Ne(840)  problem  in  the  spectral  Bent  model  described  later.  Two  tests  were 
generated  from  each  model  ionosphere  -  a  single-satellite  test  with  data  along 
a  simulated  DMSP  F6  satellite  orbit  (nominally  dawn-dusk)  and  a  two-satellite 
case  with  data  along  simulated  F6  and  F7  (nominally  10LT  -  22LT)  orbits. 
These  eight  model  tests  were  to  be  supplemented  by  daily  SSIELD  runs  on 
whatever  SSIE  data  from  satellite  F7  were  available  in  the  AGDB. 

Prior  to  the  actual  testing,  several  coding  errors  were  discovered  in 
extracting  the  spectral  Bent  model  routines  from  SSIELD,  as  well  as  severe 
problems  with  several  of  the  model  coefficient  sets  used  to  represent  the  Bent 
N  (840)  values.  The  coding  errors  were  corrected  prior  to  testing,  but 
correction  of  the  special  Bent  coefficients  were  beyond  the  scope  of  the 
assessment.  Figures  5  and  6  illustrate  this  problem  with  the  Bent 

coefficients.  Not  only  is  the  latitudinal  behavior  unrealistic,  but  the  model 
can  produce  negative  density  values  at  high  latitudes  during  almost  half  of 
the  year.  Since  this  problem  could  not  be  corrected  before  the  assessment 
testing,  all  model  tests  were  run  using  a  date  of  15  June  (day  166)  which,  as 
can  be  seen  in  Figure  6,  should  minimize  the  impact  of  this  problem.  It  must 
be  kept  in  mind,  however,  that  the  test  results  will  only  indicate  how  well 
SSIELD  performs  during  the  short  time  of  the  year  that  the  spectral  Bent  model 
is  reasonably  behaved  (roughly  15  June  ±  two  weeks)  and  will  not  indicate  how 
poorly  it  might  perform  at  other  times. 

4.1.3  Assessment  Test  Results  and  Evaluation 

Due  to  untimely  degradation  of  the  SSIE  sensors  on  the  F7  satellite,  no 
data  were  available  for  live  testing,  so  all  testing  was  on  bogus  data  as 
described  in  4.1.2.  The  overall  results  of  the  eight  tests  are  summarized  in 
Table  2. 
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Figure  5.  Latitudinal  variation  of  Me(840)  at  local  noon 
from  the  spectral  Bent  model . 
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TABLE  2.  SSIELO  Test  Analysis  Summary. 


Input  Model 

RMS 

Error 

RMS  % 
Error 

0600LT 

-  Beta  - 

1800LT  000LT 

2200LT 

Bent;  F6  only 

1.06xl04 

22.6 

F 

0.086 

Bent;  F6,  F7 

l.llxlO4 

22.6 

F 

0.086 

0.124 

F 

DHR(10)*;  F6  only 

1.45xl04  : 

105.2 

F 

0.144 

- 

- 

DHR(10) ;  F6,  F7 

1.48xl04 

88.1 

F 

0.144 

F 

F 

DHR(60) ;  F6  only 

2. 15xl04 

45.2 

0.086 

F 

- 

- 

DHR(60) ;  F6,  F7 

2. 12xl04 

39.4 

0.086 

F 

F 

0.141 

DHR ( 120) ;  F6  only 

3.50xl04 

33.7 

F 

F 

- 

- 

DHR( 120) ;  F6,  F7 

3.37xl04 

31.9 

F 

F 

0.237 

0.186 

The  first  two  columns  for  each 

test  case  are  the  total 

RMS  difference  and  RMS 

percent  difference 

between  the 

NJ840) 

grids 

output  by  SSIELD 

and  the  grids 

generated  by  the 

various  models  from 

which 

the  simulated  DMSP  data  were 

extracted.  The  last  four  columns  indicate  whether  the  covariance  function 
determination  (subroutine  COVEST)  was  successful  for  the  given  LT  map.  If  a 
value  is  given,  it  corresponds  to  the  B  variable  in  the  covariance  function; 
an  "F"  indicates  that  COVEST  failed.  Since  the  analyses  for  0600LT  and  1800LT 
are  the  same  for  either  one-  or  two-satellite  inputs,  the  B  values  for  these 
times  will  be  the  same  in  both  cases,  for  a  total  of  16  independent  LT  map 
analyses  in  the  8  test  cases.  Of  these  16  analyses,  the  COVEST  analysis 
failed  9  times. 

Figure  7  shows  a  breakdown  of  the  RMS  values  given  in  Table  2  as 
functions  of  latitude  and  local  time,  respectively,  for  the  4  two-satellite 
tests.  The  solid  lines  in  these  figures  are  for  the  DHR  model  inputs,  the 
dashed  lines  are  for  the  spectral  Bent  model  inputs.  The  problems  at  80° 
latitude  are  primarily  due  to  inexplicably  high  Ng(840)  values  in  the  SSIELD 

The  numbers  in  parentheses  indicate  the  sunspot  number  used  by  the  model. 
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Figure  7.  Variation  of  RMS  error  and  RMS  percent  error  with  latitude  and  with  local 
time  for  all  two-satellite  tests. 


analysis  at  IOOOlT  and  1800LT,  but  the  equatorial  increase  in  RMS  percent 
error  occurs  at  all  local  times.  The  local  time  plots  show,  rather 
graphically,  definite  problems  with  the  LT-spreadinq  algorithm.  The  vertical 
dotted  lines  are  at  local  times  where  SSIELO  constructed  LT-maps  (0600,  1000, 
1300,  and  ??00LT) .  The  fact  that  something  "happens"  in  the  analysis  at  these 
times ,  particularly  at  1800LT,  indicates  a  failure  of  the  LT-spreading 
process.  This  can  also  ne  seen  in  F  gure  8,  wnich  shows  the  longitudinal 
variation  of  N  (840)  for  three  latitudes  from  the  00  LIT  map  of  the  DHR(60)  two- 
satellite  test.  Again,  the  analysis  shows  definite  discontinuities  at  the 
times  where  LT  maps  were  constructed,  indicated  on  the  figure  by  vertical 
dotted  lines.  (Note:  The  RMS  error  summary  plot  in  Figure  S  does  not  show  a 
dip  at  IOOOLT  due  to  purely  coincidental  balancing  of  smaller  errors  at 
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Figure  8.  Representative  results  of  the  SS I  ELD  assessment  tests. 


latitudes  00°-70°  and  a  very  large  error  at  80°.  This  same  effect  caused  the 
increases  instead  of  dips  in  the  RMS  percent  error  plot  at  10001T  and  1800LT.) 

This  high-latitude  problem  is  caused  somewhere  in  the  construction  of  the 
IT  maps.  Figures  9a-d  show  the  variation  of  Ng(840)  with  latitude  at  0600, 
1800,  1000,  and  2200LT  from  the  OOUT  map  of  the  DHR(60)  two-satellite  test. 
The  solid  lines  are  from  the  model,  the  dashed  lines  are  results  from  the  LT 
maps,  and  the  dotted  lines  are  from  the  final  SSIELD  output.  The  large 
increase  at  80°  latitudes  shows  up  quite  clearly  in  the  1000  and  1800LT  maps. 
The  problem  at  the  equator  also  shows  up  in  these  plots,  indicating  that  it, 
too,  may  be  due  to  a  problem  in  the  LT  map  analysis. 
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Fiqure  9.  Latitudinal  variation  of  Ne(640)  from  the  DHR(60)  data,  the  LT  map 
analysis  output,  and  the  final  SS I  ELD  analysis  output. 


In  summary,  the  following  problems  and  deficiencies  were  identified: 

1.  LT  map  construction. 

a.  The  decimation  scneme  used  to  reduce  the  input  data  to  100 
points  per  LT  map  effectively  throws  out  over  60%  of  tne  availaole  northern- 
hemisphere  data. 

b.  Interpolation  in  time  is  not  employed  in  using  the  spectral 
Bent  model  to  "correct"  data  to  the  desired  local  time,  (see  3b.  below). 


c.  The  assumptions  made  in  constructing  the  covariance  function 
that  the  field  is  homogeneous  and  isotropic  may  be  unwarranted.  In 
particular,  studies  have  found  that  correlation  distances  of  ionospheric 

density  variations  are  different  in  the  magnetic  north-south  and  east-west 
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directions  . 

d.  The  B -parameter  estimation  scheme  failed  to  produce  a  solution 
in  9  of  16  test  cases. 

e.  The  LT  maps  fit  the  data  poorly  at  the  "edges"  of  the  field, 

i.e.,  at  high  latitudes  and  near  the  equator. 

2.  Time  interpolation. 

a.  The  model  used  for  interpolating  can  produce  unrealistic 
variations  in  Ng ( 840 ) . 

b.  The  interpolation  scheme  can  produce  very  sharp 
discontinuities  at  times  when  data  (from  the  LT  maps)  are  available, 
particularly  if  the  data  values  are  far  different  from  those  calculated  by  the 
model  used  for  interpolating. 

c.  If  the  LT  map  analysis  fails  to  provide  an  Ng(840)  value  for  a 
particular  grid  point,  a  value  of  zero  is  generated  for  that  grid  point  rather 
than  some  representative  Ng(840)  value. 

3.  Spectral  Bent  model. 

a.  Fields  of  Ng(840)  produced  by  the  model  can  be  unrealistic, 
occassionally  containing  negative  density  values. 

b.  The  model  is  discontinuous  in  time.  A  different  set  of 
coefficients  are  used  for  each  of  24  local  times  with  no  interpolation  for 
times  not  on  the  hour. 

c.  The  model  is  valid  only  for  a  single,  unstated,  solar  epoch.  A 
new  model  would  need  to  be  implemented  as  the  solar  epoch,  as  indicated  by  the 
mean  sunspot  number,  changes. 


4.  General 


a.  The  program  has  no  "memory"  capabi 1 ity;  i .e. ,  the  starting  point 
of  the  analysis  is  always  the  internal  model. 

b.  Expanding  the  program  to  provide  other  capabilities,  such  as 
including  the  southern  hemisphere  in  the  analysis,  would  be  difficult. 

c.  The  documentation  provided  requires  additional  work  to  meet 
AFGWC  and  DoD  standards,  and  neither  the  program  code  nor  the  in-line 
documentation  (comments)  meets  AFGWC  standards. 

The  limitations  and  shortcomings  of  the  current  SSI  ELD  program  were 

discussed  with  Dr.  Rich  (AFGl/PHG)  and  Mr.  Bussey  ( AFGWC/TSIS) ,  and  a 

consensus  was  reached  that  although  the  analysis  methods  used  in  SS I  ELD  were 
theoretically  sound,  they  were  probably  inappropriate  for  the  current 
application  and  would  be  inadequate  to  tr,c  task  even  if  the  problems  were 
corrected.  It  was  decided,  therefore,  to  redirect  the  near-term  efforts  under 
tnis  contract  to  identifying  and  investigating  alternative  analysis 
techniques  for  producing  gridded  fields  of  SSIE  Ng(840)  data. 

Toward  this  end,  informal  discussions  were  held  with  Dr.  Rich  and  Mr. 
Bussey  to  determine  the  requirements  be  met  by  the  SSIELD  program.  A 

distillation  of  these  discussions  provided  the  following  list  of  requirements 

and  desired  attributes: 

1.  The  program  must  be  documented  and  tested  in  accordance  with  DoD 
Standard  7935,  dated  15  Feb  1983?5,  and  with  the  AFGWC  software  standards^. 

2.  The  program  must  be  coaed  in  ANSI  Standard  X3.9-1978,  as 
implemented  in  the  UNI  VAC  ASCII  FORTRAN. 

3.  Any  internal  model  of  Nc(840)  usea  should  oe  a  continuous  function  of 
latitude,  longitude,  time,  season  ( day-of -year ) ,  and  sunspot  number  (solar 
epoch) . 

4.  Ideally,  any  internal  N  (84C)  moc’l  should  oe  oasec  on  satellite 
observations.  Provisions  for  updating  lv.s  model  should  be  included. 
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5.  Data  from  the  SSIE  sensor  should  be  fetched  from  a  user-defined 
section  of  the  Astrogeophysical  Data  Base  (AGDB). 

6.  Provisions  for  future  incorporation  of  TEC  data,  from  either  the 
AGDB  or  from  a  TEC  data  preprocessor,  should  be  included. 

7.  The  analysis  algorithm  should  use  all  data  available  and  be  able  to 
process  data  from  both  northern  and  southern  hemispheres. 

8.  If  no  data  are  available  for  a  particular  program  run,  the  analysis 
should  default  to  the  internal  model. 

9.  The  analysis  should  incorporate  a  "memory"  capability;  i.e.,  the 
analysis  field  from  a  previous  run  could  be  used  as  the  starting  point  for  the 
current  run. 

10.  The  final  analysis  field  should  not  contain  distinct  features 
(bull's-eyes  or  ridges)  at  points  where  data  are  available. 

11.  The  locations  of  the  grid  points  for  the  output  grids  should  be 
defined  such  that  they  may  be  easily  changed. 

Although  this  is  not  an  officially  approved  statement  of  requirements,  it 
provides  a  starting  point  for  an  investigation  of  alternative  analysis 
methods  for  incorporation  in  program  SSIELD. 

4.2  Task  4:  Ne(840)  Data  Preprocessor  Improvements 

The  starting  point  for  this  study  was  a  general  review  of  various 

analysis  methods  developed  for  use  in  other  ionospheric  applications,  in 
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particular,  those  developed  in  constructing  the  ITS78  model  ionosphere  ’  , 

those  used  in  AFGWC  program  UKFIIE  for  updating  the  ITS78  coefficients^,  and 
the  global  analysis  method  used  in  the  4D  model^.  From  this  review,  and  the 
requirements  listed  in  Section  4.1.3,  the  following  general  design  was 
developed  for  an  improved  SSIELD  program: 

1.  The  core  of  the  program  would  be  a  spectral  model  of  Ne(840)  which 
would  include  continuous  variability  with  location,  time,  season,  and  solar 

epoch.  The  format  of  this  model  would  be  similar  to  the  ITS78  f  F2  and  M3000 

13  0 

spectral  models  . 

2.  Input  SSIE  Ng(840)  data  would  be  used  to  modify  the  coefficients  of 
this  model  in  two  stages: 


a.  All  mid-latitude  data  would  be  used  to  adjust  the  model 

through  derivation  of  an  effective  10.7cm  solar  flux  (FIG)  defined  as  the  F10 

value  which,  when  used  by  the  Ne(840)  model,  minimizes  the  mean  difference 

between  the  model  and  the  observations.  A  similar  method  for  adjusting  the 
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ITS78  fQF2  coefficients  has  been  used  at  AFGwC  for  several  years  . 

b.  All  data  would  be  used  to  modify  tne  coefficients  resulting 
from  the  effective  F10  analysis,  or  from  a  previous  analysis,  using  a  spectral 

?Q 

data  assimilation  method  developed  by  Flattery"  and  used  in  the  global 
analysis  section  of  the  4D  model7*. 

3.  The  required  UT  grids  would  be  generated  from  either  the  base 
Ne(840)  model,  the  results  of  the  effective  F10  analysis,  the  results  of  the 
global  analysis,  or  from  the  results  of  some  previous  analysis,  at  the 
discretion  of  the  user. 

Investigation  of  this  design  would  be  accomplished  in  four  phases: 

1.  Define  the  form  of  the  spectral  N  (840)  model.  This  would  include 
selection  of  coordinate  system,  functions  to  be  used  in  the  spectral 
expansion,  truncation  limits  on  the  spectral  expansion,  and  selection  of  an 
initial  model  of  Ne(840)  on  which  to  base  the  spectral  model.  The  final  step 
of  this  phase  would  be  the  generation  of  coefficients  for  the  spectral  Ng(840) 
model . 

?.  Develop  the  effective  F10  analysis  technique.  The  starting  point 
for  this  phase  would  be  the  AFGWC  UKFILE  program^. 

3.  Develop  the  global  assimilation  technique.  The  starting  point  for 
this  phase  would  be  the  global  analysis  section  of  the  4D  model'7. 

4.  Test  the  full  analysis  program  on  both  model  and  "real-world" 
inputs  to  determine  how  best  to  control  the  analysis  so  the  occurance  of 
spurious  features  at  data  points  is  reduced  with  minimum  loss  in  analysis 
accuracy . 

work  nas  begun  on  phases  one  and  t.v-  or  mis  study. 
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4.2.1  Ne(840)  Spectral  Model  Definition 

The  first  step  in  defining  the  spectral  N0(84O)  model  was  to  select 

coordinate  systems  in  which  to  describe  the  model's  variations  in  solar  epoch, 

season,  time-of-day,  and  location.  It  was  decided  to  use  the  same  coordinates 

as  the  ITS78  model  for  solar  epoch  (F10  -  the  solar  10.7cm  radio  flux  for 

season  (DOY  -  the  day-of-the-year) ,  and  for  time-of-day  (UT  -  universal 

time)),  but  to  use  a  different  system  for  location.  In  developing  the  ITS78 
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model,  Jones  et  ai  found  that  using  a  coordinate  system  based  on  the 
configuration  of  the  earth's  magnetic  field  improved  the  analysis.  Their 
choice  of  location  coordinates  was  the  modified  dip  latitude  and  geographic 
longitude,  the  modified  dip  latitude  being  defined  as 


x  =  tan 


/cos  X 


where  I  is  the  magnetic  dip  angle  and  X  is  the  geographic  latitude.  This 
system  was  not  selected  for  this  application  as  it  did  not  meet  the  following 
requirements: 

1.  The  system  should  provide  both  magnetic  latitude  and  longitude 
within  the  same  conceptual  framework,  as  well  as  a  logical  extension  to  the 
calculation  of  a  magnetic  local  time. 

2.  The  latitude  coordinate  should  behave  as  dip  latitude  near  the  dip 
equator  and  as  corrected  geomagnetic  (CG)  latitude  in  the  auroral  and  polar 
regions. 

3.  Both  geographic-to-magnetic  and  magnetic-to-geographic  transfor¬ 
mations  are  needed,  and  should  be  computationally  fast. 

An  extensive  study  of  magnetic  coordinates  and  their  applicability  to 

ionospheric  mapping  problems  conducted  under  another  project  led  to  the 
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development  of  F-layer  apex  coordinates  ,  a  modification  of  apex  coordinates 
developed  by  Van  Zandt  et  ai  .  Figure  10  is  a  contour  map  of  F-layer  apex 
latitude  and  logitude.  The  latitude  spacing  is  10°,  and  the  longitude  spacing 
is  90°.  The  horizontal  and  vertical  dashed  lines  indicate  the  F-layer  apex 
equator  and  zero  meridian  respectively.  This  system  is  defined  similarly  to 
CG  coordinates,  but  with  respect  to  the  dip  equator  instead  of  the  dipole 
equator.  Thus,  since  the  dip  and  dipole  equators  merge  at  large  distances 
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Figure  10.  F-layer  apex  coordinates. 


from  the  earth,  F-layer  apex  latitude  is  similar  to  dip  latitude  near  the 

common  F-layer/dip  equator  and  becomes  identical  to  CG  latitude  at  high 

latitudes.  Algorithms  to  provide  the  required  coordinate  transforms  were 

developed  previously,  as  well  as  a  look-up  table  for  rapid  calculation  of 
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r-layer  apex  coordinates  from  geographic  coordinates  .  Since  this  system 
meets  all  the  stated  requ i rements ,  and  is  readily  available,  it  was  decided  to 
use  F-layer  apex  coordinates  for  the  location  coordinate  system. 

The  next  step  in  the  model  definition  was  to  select  the  expansion 
functions  for  each  model  variable.  As  with  the  selection  of  the  coordinate 
systems,  it  was  decided  to  use  the  functions  selected  for  the  ITS78  model  for 
the  solar  epoch,  season,  and  time-of-day  variations.  The  solar  epoch 
variation  is  expressed  as  a  polynomial  function  in  FiO,  tne  season  variation 
as  a  trignometrio  series  in  DGY,  unc  tna  time-of -day  variation  as  a 

97 

trignometnc  series  hi  U i  .  For  ideation  functions,  Jones  et  al  used  a  set 
of  trignometr ic  polynomials  which  are  orthonormal ized  with  respect  to  the 
locations  of  tne  input  oat  a  for  each  analysis  run.  In  the  40  model, 
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orthonormal  associated  Legendre  functions,  which  are  very  similar  to  the  ad- 
hoc  set  used  in  ITS78,  are  used  to  model  the  latitude-longitude  variation.  It 
was  decided  to  use  the  associated  Legendre  functions  for  the  following 
reasons: 

1.  The  use  of  associated  Legendre  functions  for  spherical  harmonic 
analysis,  is  a  standard,  we  1 1 -documented  analysis  technique. 

2.  Associated  Legendre  functions  were  chosen  for  use  in  the  4D  model  to 
avoid  singularity  problems  encountered  in  using  the  ITS78  functions  in 
program  UKFILE,  the  predecessor  the  4D  model^. 

3.  The  method  to  be  used  for  the  global  assimilation,  the  4D  model 
global  analysis  method,  is  designed  around  a  field  expressed  in  associated 
Legendre  functions.  Thus,  the  amount  of  design  work  required  is  reduced  if 
associated  Legendre  functions  are  used. 

The  final  form  selected  for  the  Ng(840)  model  is  as  follows: 

L  M  N 

NeO  ,  v,t,F  ,D;840)  =  £  2  2m  Yffmn(F,D)  T  (t)T(*)  P*/2(sin\)  (9) 

£  =  1  m=l  n=  — 


mn 


(F,D)  =  £ 
r=0 


2 

i=0 


mn. 


r,T,(0) 


(10; 


where  A  and  $  are  F-layer  apex  latitude  and  longitude;  t 
solar  flux;  D  is  the  day-of-the-year ;  P^^(sinA) 
associated  Legendre  functions  of  degree  n  and  order 
T.j(D)  are  orthonormal  trignometric  series  given  by 


is  UT;  F  is  the  10.7cm 
are  the  orthonormal 
and  T  (t),  TJ  :),  and 


T2j-i(x>  =  — =—  cos  jx  , 

T  0  •  ( x )  =  -4: 


sin  jx 


The  next  step  is  to  select  values  for  the  five  truncation  limits  in  equations 
(9)  and  (10),  L,  M,  N,  R,  and  I. 


The  selection  of  the  five  truncation  limits  involves  a  rather  com¬ 
plicated  trade-off  analysis  between  analysis  accuracy  and  resolution  require¬ 
ments  and  operational  size  and  run-time  constraints.  As  no  statement  of 
operation  constraints  is  available,  the  following  targets  were  established 
based  on  informal  discussions  witn  AFGWC  personnel:  1)  run-time  on  the  oraer 
of  four  minutes  of  CPU  time  or  less  on  a  UN  I  VAC  1100/82  for  a  24-hour,  two- 
satellite  analysis,  2)  executable  program  size  under  6GK,  and  3)  disk  space 
requirements  on  the  order  of  one  position  (  -  115K  words)  or  less  for  all 
required  data  files.  It  was  decided  to  approach  the  accuracy/resolution 
question  by  determining  the  maximum  resolution  afforded  by  the  operational 
constraints.  First,  the  truncation  limits  for  equation  (10)  were  set  as  the 
same  as  those  used  in  the  ITS78  model  -  1=8  and  R=2,  as  was  the  value  of  the 
truncation  limit  for  the  UT  expansion  in  equation  (9)  -  1=13.  This  leaves  the 
truncation  limits  on  the  associated  Legendre  function  expansion,  N,  and  M. 

There  are  several  truncation  schemes  available  for  a  finite  associated 
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Legendre  function  series  ,  such  as  the  rhomboid  scheme  used  in  the  4D  model, 
shown  graphically  in  Figure  11.  A  modified  triangular  scheme,  shown  in  Figure 
12,  was  selected  for  the  Ne(840)  model  based  on  the  following  considerations: 

1.  The  truncation  scheme  inherent  in  the  ITS78  model,  also  shown  in 
Figure  12,  approximates  a  modified  triangular  truncation. 


2.  This  scheme  requires  fewer  functions  than  the  4D  model  scheme  for 
the  same  maximum  function  degree  for  low  longitudinal  wave  numbers. 

3.  This  scheme  acts  as  a  filter  to  smooth  out  features  which  are  small 
in  both  latitudinal  and  longitudinal  extent. 


The  limits  used  in  the  truncation,  M  =10*  and  N=23  at 


in  =0*, 
P 


were  selected  to 


*The  truncation  1  imit,  M  ,  and  index,  m  useu  in  this  discussion  are  related  to  M 

r  P 

and  m  in  equation  ( 9)  by 
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Figure  11.  The  rhomboid  truncation  scheme  for  an  associated  Legendre  function 
series  used  in  the  4D  model. 
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Figure  12.  The  modified  triangle  truncation  scheme  (heavy  line)  chosen  for  the 
Ne(840)  model  and  the  inherent  truncation  employed  in  the  ITS 7S  model 
(fine  line). 


match  at  least  the  value  for  N  at  m  =0  used  in  the  4D  model,  N=19.  The  limits 
in  equation  (9)  for  this  truncation  are  M=21  and  N=20  - 

A  preliminary  design  for  the  data  file  to  be  used  by  this  program  was 
developed  to  estimate  the  disk  space  it  would  require.  Included  in  the  file 
would  be  the  F-layer  apex-coordinate  transformation  table  (8190  words,  293 
28-word  sectors*),  the  coefficients  for  the  spectral  Ng(840)  model  (one  sector 
per  time-location  coefficient,  2873  sectors  total),  five  sets  of  coefficients 
from  previous  analyses  (2373  words,  103  sectors),  and  six  file-control 
sectors,  for  a  total  of  3687  sectors.  This  is  roughly  103K  words,  well  within 
the  target  of  115K. 

The  final  step  in  defining  the  spectral  Ne(840)  model,  prior  to 
generating  the  inoue'l  coefficients,  was  to  select  a  model  of  N  (840)  from  which 
to  generate  the  coefficients.  Since  this  analysis  program  is  being  developed 
for  use  at  AFGWC,  it  was  decided  to  use  the  DHR  profile  model  in  conjunction 

with  the  ITS78  f  F2  and  M3000  models,  at  least  for  initial  design  and  testing. 
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It  is  planned  to  obtain  a  copy  of  the  Bent  profile  model  and  construct  a  set 

of  spectral  model  coefficients  for  it  as  well,  for  comparison  testing  with 

actual  SSIE  data.  The  software  used  to  calculate  the  spectral  model 

coefficients  and  all  routines  which  use  them  will  be  designed  such  that  they 

are  completely  independent  of  the  Ng(840)  model  used  in  generating  the 

coefficients. 

The  coefficients  of  the  spectral  N  (840)  model  are  calculated  by  a  least- 
squares  fit  of  the  model  functions  to  fields  generated  by  the  ITS78-DHR  (or 
ITS78-Bent)  models.  In  order  to  reduce  the  dynamic  range  of  the  Ng(840) 
fields,  the  probable  source  of  the  negative  Ng(840)  values  in  the  SSIELD 
model,  the  density  values  are  converted  to  the  equivalent  plasma  frequency, 
f  (840),  prior  to  performing  the  least-squares  fits.  The  ITS78-DHR  model  is 
used  to  build  10  ’  latitude  by  15°  longitude  grids  in  F-layer  apex  coordinates 
covering  the  range  -80°  to  +80°  latitude  for  each  UT  hour  (00-23),  for  each 
month  of  the  year,  for  F10  values  of  70,  120,  170,  and  220.  Each  of  the  48 
sets  of  time-location  grids  is  reduced  to  48  sets  of  2873  coefficients  using 
least-squares  fits  to  the  time,  T.(t),  and  location,  T( ; ) Pn( sin  A ) , 
functions  in  equation  (9).  These  are  further  reduced  to  12  sets  of  8619 

*The~  FTTe'sTze  iV  defined  in  terms  of  28  word  sectors  as  the  UN  I  VAC  1100/82 
computer  uses  sector-addressable  1/0  to  disk  mass  storage.  Disk  space  is 
allocated  by  track  (64  sectors)  or  position  (64  tracks). 


coefficients  by  a  least-squares  fit  in  F10,  and  finally  to  a  set  of  77,571 
coefficients  by  a  fit  in  season. 

A  complete  coefficient  set  has  not  yet  been  generated,  although  sub-sets 
have  been  generated  to  test  all  three  fitting  procedures.  Figures  13  and  14 
show  sample  results  of  the  time-location  fits.  The  upper  plot  in  each  figure 
is  the  ITS78-DHR  Ne(840)  field,  and  the  lower  plot  is  from  equation  (9)  using 
the  coefficients  from  the  time-location  fit.  These  plots  are  in  F-layer  apex 
coordinates  with  a  geographic  latitude-longitude  grid  overlay.  Figures  15 
and  16  are  examples  of  the  results  of  the  F10  fits,  showing  sample  latitude 
and  UT  variations,  and  Figure  17  shows  the  seasonal  fit  to  the  first  time- 
location  coefficient  for  F10  =  70.  All  of  the  fits  are  within  acceptable 
limits.  A  complete  set  of  coefficients  modeling  the  F10  and  time-location 
variations  for  the  month  of  May  Have  been  generated  for  initial  testing  of  the 
effective  F10  analysis. 


4.2.2  Effective  F10  Analysis 

The  first  stage  in  the  new  analysis  scheme  consists  of  adjusting  the 
value  of  F10  used  in  the  coefficient  equation 

y  r  -  E  s  T  f ior 


where 

9 

emn  =  £  a*mn 
r  t—4  ri  i  ' 
i  =  l 

until  the  average  difference  between  all  mid-latitude  Ne( 840)  data  and  the 

values  calculated  using  equation  (9)  is  minimized.  This  value  of  F10  is 

denoted  the  effective  F10  value  for  the  analysis  and  is  conceptually 

equivalent  to  the  effective  sunspot  number  (SSN)  calculated  in  AFGWC  program 
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UKFILE  by  fitting  observed  fQF2  values  with  the  ITS78  model  .  An  initial 
test  of  this  technique,  using  the  binary  search  algorithm  from  UKFILE,  was  run 
using  data  along  a  simulated  dawn-dusk  DMSP  orbit  extracted  from  a  time- 
location  Ne(840)  grid  set  calculated  for  F10  =  170.  Given  an  initial  F10 
value  of  120,  the  analysis  calculated  an  effective  F10  of  170.2  in  6 
iterations.  Althouoh  the  accuracy  shown  in  this  test  was  more  than 

acceptable,  it  is  hoped  that  using  a  more  sophisticated  search  algorithm  could 
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Figure  17.  Sample  results  from  the  season  analysis.  Shown  is  a  plot  of 
the  Yq°  coefficient  (see  Equation  10)  from  each  month  of  the 
time- location  analysis  for  F10-70  (crosses)  and  the  results 
_ of  the  season  analysis  for  this  coefficient. 


possibly  achieve  similar  accuracy  in  fewer  iterations.  An  evaluation  of  other 
search  algorithms  is  the  next  step  in  the  current  study.  When  this  is 
completed,  work  will  begin  on  the  global  analysis  method. 


4.3  Task  5:  Combined  TEC  and  Ng(840)  Data  Preprocessor 

Other  than  informal  discussions  with  personnel  from  AFGWC/TSIS  (Mr.  Bob 
Bussey)  and  AFGL/LIS  (Mr.  Jack  Klobuchar),  little  has  been  done  on  this  Task 
pending  resolution  of  the  problems  with  the  current  Ng(840)  preprocessor  and 
formal  tasking  within  AFGWC  for  Mr.  Bussey  to  work  on  improvements  to  the  TEC 
preprocessor.  One  possible  method  for  combining  the  two  data  types  in  a 
single  analysis  would  be  to  convert  TEC  to  an  effective  Ne(840)  and  process  it 
with  the  SSIE  Ng(840)  data  .  This  will  be  investigated  when  the  work 
described  in  Section  4.2  has  been  completed. 
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5.  4D  MODEL  IMPROVEMENTS  STUDY 

The  two  tasks  which  make  up  this  study,  Tasks  1  and  3,  were  both 
completed  during  this  project  year,  although  Task  3  was  modified  slightly  as  a 
result  of  the  IRI79  model  study  discussed  earlier. 

5.1  Task  1:  4D  Model  Global  Analysis 

The  4D  model,  as  originally  configured,  produced  analyses  of  the 
northern  hemisphere  only  in  order  to  reduce  both  program  size  and  model  run¬ 
time.  This  is  done  by  assuming  that  the  southern  hemisphere  is  a  mirror-image 
of  the  northern  hemisphere,  thus  forcing  the  coefficients  of  all  odd  asso¬ 
ciated  Legendre  functions,  P1^,  to  be  zero.  Thus,  only  the  even  P™  functions 
and  their  coefficients  need  be  processed,  reducing  the  total  number  of 
coefficients  by  a  factor  of  two,  as  well  as  the  run  time  of  the  global  analysis 
section  of  the  4D  model.  Unfortunately,  the  result  was  that  half  of  the  data 
collected  by  the  SSIE  sensor,  i.e.  all  southern  hemisphere  data,  were  not  usable 
by  the  4D  model.  The  objective  of  Task  1  was  to  determine  what  modifications 
were  required  to  configure  the  4D  model  for  a  true  global  analysis. 

As  it  turned  out,  only  a  few  modifications  to  the  4D  model  are  required 
to  allow  the  model  to  be  configured  for  either  global  or  hemispheric  analysis. 
These  were: 

1.  Insure  that  the  routine  used  to  calculate  the  associated  Legendre 
functions  calculates  both  even  and  odd  functions. 

2.  Set  up  a  parameter  which  specifies  the  model  configuration. 

3.  Set  up  parameters  for  all  array  dimensions,  array  sizes,  and  DO- 
LOOP  range  variables  to  be  functions  of  the  configuration  parameter. 

4.  Change  all  code  in  the  model  which  limits  the  analysis  to  the 
northern  hemisphere  to  be  configuration  dependent. 


5.  Document  the  procedures  for  configuring  the  4D  model  in  full  global  !] 

analysis  mode.  I 


I 
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These  changes  were  implemented  in  the  40  model  as  part  of  a  separate  contract, 
and  the  configuration  change  procedures  were  included  in  the  documentation 

g 

provided  as  part  of  that  project  . 

These  modifications  will  improve  the  utilization  of  southern  hemisphere 
SSIE  observations,  but  it  is  doubtful  that  this  mode  will  be  operationally 
useable  for  real-time  applications,  as  it  doubles  both  the  number  of  coeffi¬ 
cients  required  for  analysis  and  the  model  run-time  over  a  single-hemisphere 
analysis. 

5.2  Task  3:  Model-based  4D  Height  Functions 

The  4D  model  uses  a  set  of  empirically-derived  orthonormal  functions  to 

describe  the  height  variation  of  the  ionospheric  electron  density'7.  The 

original  function  set  was  derived  from  a  set  of  electron  density  profiles 

collected  by  the  Millstone  Hill  incoherent  scatter  radar,  extrapolated  down- 

18 

ward  to  95km  and  upward  to  2000km  using  the  OHR  model  .  Limited  studies  of 

IQ  ?n 

the  4D  model  during  the  late  1970's  ’  showed  that  the  Millstone  Hill 

function  set  could  not  adequately  match  observed  electron  density  profiles  at 
other  than  mid-latitudes,  and  for  moderate  to  high  sunspot  numbers  at  all 
latitudes.  This  problem  was  magnified  by  the  differences  in  shape  between 
profiles  built  by  the  DHR  model  from  i ..put  data  and  those  produced  by  the 
Millstone  Hill  functions. 

The  optimal  set  of  data  from  which  to  construct  a  more  representative 
height  function  set  would  be  a  collection  of  incoherent  radar  soundings 
covering  the  height  range  l00-2000km  from  equatorial,  middle,  auroral,  and 

polar  latitudes  at  various  longitudes  for  at  least  one  solar  cycle.  As  the 

collection  of  such  a  data  set  is  unfeasible,  it  was  proposed  to  use  the  IRI79 
model  to  construct  such  a  data  set  from  which  height  functions  could  be 
constructed.  The  reasons  for  choosing  the  IRI79  model  were  twofold: 

1.  The  IRI79  model  was  constructed  from  as  extensive  a  data  base  as 
possible  to  cover  most  ionospheric  regimes,  and 

2.  It  was  proposed  to  incorporate  the  IRI79  profile  parametrization 
into  the  40  model  (Task  2  of  this  project).  Thus,  there  would  be  no  conflicts 
between  the  internal  profile  model  and  the  height  functions. 


The  results  of  the  IR 179  study  conducted  in  Task  2,  and  the  decision 
based  on  that  study  not  to  use  the  I R I 79  parametrization  within  the  4D  model, 
have  led  to  a  modification  of  this  task.  Since  the  4D  will  continue  to  use  the 
DHR  profile  model,  it  was  decided  to  construct  height  function  sets  from  a 
collection  of  profiles  built  from  the  DHR  model.  While  this  is  a  somewhat 
less  than  optimal  data  set,  this  will  at  least  alleviate  problems  caused  by 
conflicts  between  the  internal  profile  model  and  the  height  functions. 

Three  sets  of  DHR  model  profiles  were  generated  and  height  function  sets 
calculated  from  them  during  visits  to  AFGWC.  All  three  sets  consisted  of  a 
collection  of  5184  profiles  generated  for  a  10°  latitude  by  15°  longitude  by 
lh  local  time,  northern-hemisphere  grid  covering  a  height  range  of  95-2000km. 
Set  1  included  profiles  generated  for  four  solar  epochs  (sunspot  numbers  10, 
50,  90,  and  130)  and  two  seasons  (15  June  and  15  December),  set  2  included 
profiles  generated  for  two  solar  epochs  (sunspot  numbers  20  and  120)  and  four 
seasons  (15  March,  15  June,  15  September,  and  15  December),  and  set  3  was 
generated  for  a  single  sunspot  number  (56)  and  date  (20  September).  Figures 
18a-18d  show  the  four  height  functions  extracted  from  these  three  sets,  as 
well  as  the  original  Millstone  Hill  functions.  Only  the  lower  1000  kilometers 
are  plotted,  as  all  four  functions  shown  decay  exponentially  toward  zero  in 
the  1000-2000km  height  range. 

These  new  height  function  sets  and  the  software  reguired  to  generate 
additional  sets  were  provided  to  AFGWC,  and  the  procedures  for  height 
function  set  generation  were  formally  documented  as  part  of  a  separate 

g 

contract  .  Additional  sets  for  other  solar  epoch/season  combinations  will  be 
generated  during  the  next  project  year,  as  well  as  for  the  new  height  profile 
model  being  developed  under  Task  6,  should  it  prove  to  be  a  viable  replacement 
for  the  DHR  model . 
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Figure  18a.  First  40  height  functions  derived  from  Millstone 
and  three  DHR  model  profile  data  sets. 
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Figure  18b.  Second  40  height  functions  as  in  Figure  18a. 
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6.  CONCLUSION 


Work  on  this  project  has  progressed  satisfactorily  during  its  first 
year.  The  only  difficulties  encountered  were  the  unexpected  inadequacy  of  the 
I R 1 79  topside  parametrization  and  the  problems  described  in  Section  4.1  with 
the  current  SSIE  Ng(840)  preprocessor,  program  SSIELD.  The  main  efforts 
planned  for  the  second  project-year  are  continuation  of  work  on  Task  4  to 
improve  the  Ng(840)  preprocessor ,  and  on  Task  6,  to  investigate  alternative 
topside  parametr izations  for  the  40  profile  model.  At  the  completion  of  Task 
4,  work  will  begin  on  Task  5,  an  investigation  of  methods  for  combining  TEC 
and  Ne(840)  data  into  a  single  analysis. 
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